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Abstract. We analyze the dynamics of occupation probabilities for a certain type 
of design models by the use of two different methods. On the one hand we present 
some numerical calculations for two concrete interactions which point out that the 
occurrence of statistical dynamics depends on the interaction structure. Further- 
more we show an analytical derivation for an infinite system that yields statistical 
behaviour for the average over the whole ensemble of interactions in the Van 
Hove- limit. 



1 Introduction 



The emergence of statistical behaviour from microscopic dynamics is of special interest for two 
reasons. Its existence is evident from countless experiments, nevertheless its explanation seems 
subtle. And if statistical dynamics are established, their description is much simpler than the 
description of microscopic dynamics. From first principles the dynamics of quantum systems 
are controlled by the Schrodinger equation. Nevertheless, it has been observed that statistical 
relaxation may appear in such systems for certain quantities under certain conditions. The 
dynamics of the quantities P„ are called statistical if they are given by a master equation of 
the form 



This is a set of coupled rate equations whose solutions decay exponentially in time. 

In this contribution the dynamics for a certain type of modular design model are approached 
from two different sides. On the one hand, we show numerical calculations of the time evolution 
of occupation probabilities for a finite size version of the model with two concrete types of 
interactions. Furthermore we will present an analytical derivation that analyzes the dynamics 
of these variables for an average over all possible interactions in the Van Hove-limit and in the 
limit of an infinitely large system size, in the sense that the number of energy levels in the 
bands of the modules N goes to oo, whereas the number of modules M is always 2, i.e., M is 
still finite. 

The numerics are consistent with the analytical results for an interaction which is kept very 
general and represents the majority of all interactions. Nevertheless, there can be completely 
exceptional interactions that yield substantially different dynamics. 
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2 Numerical calculations 
2.1 Design model 

The model we investigate for the occurrence of statistical dynamics is a simple design model 
depicted in Fig. [T] (see also |ll2l3l4p . 
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Fig. 1. Sketch of the design model 



It consists of M one-dimensionally arranged identical modules that are coupled to the two 
neighbouring modules via the interaction V. Each subunit consists of a single ground state and 
an energy band. These are separated by the band gap AE. The model can be characterized by 
some rough system parameters; N is the number of levels in each band, 5e denotes the band 
width {Se AE) and A gives the average interaction strength. The Hamiltonian can be split 
into a local part Hq and an interaction V 

H = Ho + XV . (2) 

The interaction is chosen to be weak such that one can assume local energy conservation and 
restrict to the subspace of the states where one system is in the excited state and all other 
systems are in the ground state. In this contribution only chains with two modules (Af — 2) 
are regarded. V shall have an off-diagonal block form in the eigenbasis representation of Ha in 
the relevant subspace, thus the Hamiltonian is given by 
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(3) 



The relevant space is divided into two subspaces each containing the states where the /.t-th 
subsystem is in the excited state = \,2). V is normalized by 



Tr{F^} 
27V2 



= 1 . 



(4) 



One can now define the occupation probability in the /i-th subsystem P^. According to the 
predictions of the Hilbert Space Average Method (HAM) in second order (see contribution 
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Breuer/Gemmer) the time evolution of the is controUed by coupled rate equations which 
means that the dynamics are statistical 

at 

^ = RP,- RP^ , (5) 
at 

N 

with R ~ 2TrXg and 5 = ^ ■ (6) 

The rate is completely determined by the system parameters and is the same rate that appears 
in the context of Fermi's Golden Rule. In fact, HAM uses some kind of stepwise iteration of 
Fermi's Golden Rule. One also obtains some necessary conditions for the system parameters, 

Ki = X'^ « 1 , (7) 

K2 = 2A^ > 1 . (8) 

oe 

A more detailed derivation can be found in [l|2|3|4j . 

In particular, these results do not depend on the concrete realization of V. So, in principle, 
they should be valid for any possible interaction structure. 



2.2 Numerical calculations for two concrete interactions 

In this paragraph we are going to compare the dynamics of the P^ obtained from the rate equa- 
tions ([5]) with the numerical solution of the Schrodinger equation for two specific interactions. 
The system parameters are adjusted in a way that the criteria ([7]) and (H)) are well fulfilled. 

For the first interaction the matrix elements of V in the eigenbasis representation of Hq are 
chosen as random Gaussian distributed complex numbers. For the second interaction all matrix 
elements of V are set to be equal, say, Vij = 1. The random interaction possesses no structure 
at all, whereas the opposite holds true for the constant interaction. 

Figure [2] and Fig. [3] show the time evolution of Pi{t) for both interactions. The initial state 
is selected to have Pi(0) = 1. 




Fig. 2. Dynamics of Pi (t) for a completely ran- Fig. 3. Dynamics of Pi{t) for a constant in- 
dom interaction, parameters: N = 500, 5e = 0.5, teraction with Vij = 1, parameters: A'' = 500, 
A = 0.0005. Se = 0.5, A = 0.0005. 



One finds a very good agreement between both curves for the random interaction. Therefore 
the rate equations predicted by HAA/I are valid and the dynamics are statistical. For the constant 
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interaction there is no correspondence at all, so one finds no statistical relaxation in this case, 
although the system parameters, which are the same for both interactions, fulfill the criteria 
([7]) and ([8]). Referring to Sec. [3] the random interaction represents the average over all possible 
interactions very well, whereas the constant interaction demonstrates that there can be complete 
exceptions from the average. It becomes obvious that the concrete structure of the interaction 
can be indeed significant for the occurrence of statistical relaxation. 

The HAM calculation, as well as Fermi's Golden Rule, relies on second order perturba- 
tion theory, i.e., the time evolution operator is expressed by the Dyson series in second order 
truncation 

|V(r)) = (i - Jf/i(r) - ^U2{r) + ...)\m) (9) 
with the first and second order terms 

C/i(r)= r dr'Vir') , (10) 
Jo 

U2{t) = r dr' r dT"ViT')ViT") . (11) 

Jo Jo 

If this expansion is justified, all higher orders, including U2{t), must be small compared to Ui{t) 
in the relevant time regime. We use Tr{UU^} to measure the size of the respective contributions. 
HAM produces a "best guess" for the time evolution of the for a short time step (P^(i + t)) 
on the basis of Pp,{t) by using the appropriate truncation of the Dyson series, 

Pf.it + r)-P^{t)^f{T){P,,^,{t)+P,,+,{t)~2P,M) , (12) 

where fir) corresponds to a double time integral over the autocorrelation function of V (see 
also 0). One usually expects that this autocorrelation function has decayed completely after 
some decay time Tc- After Tc /(r) grows linearly with r and equation (fT2l) can be iterated 
which eventually results in the corresponding rate equations ([5]) . Figure 3] and Fig. O show a 
comparative calculation of Ui{t) and U2{t) for both analyzed interactions. Whereas Ui{t) is 
equal for both interactions, U2{t) will depend on the structure. 
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Fig. 4. Dyson terms for a completely random 
interaction, parameters: A'^ = 500, Se — 0.5, A — 
0.0005. 
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Fig. 5. Dyson terms for a constant interaction 
with Vij = 1, parameters: A'^ — 500, Se = 0.5, 
A = 0.0005. 



Ui{t) and U2{t) must be compared in the time regime where equation (|12p is iterated, i.e., 
in the region after Tc (here Tc ~ 25). This turns out to be equivalent to the beginning of the 
linear regime of Ui{t). For the completely random interaction one finds that Ui{t) ^ U2{t) at 
time Tc, so the second order truncation seems to be suggestive in this case. There is no reason 
why higher orders beyond the second should have a decisive influence, since the second order 
already leads to a very good approximation. This corresponds to the numerical calculations that 
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showed statistical relaxation. For the constant interaction U2{t) is already larger than Ui{t) 
at time Tc- So the truncation is not justified which conforms to the numerically calculated non- 
statistical dynamics. It should be remarked that this breakdown of the exponential behaviour is 
actually a problem of the structure of the interaction. It cannot be solved by simply decreasing 
the average interaction strength A. This is demonstrated in Fig. [S] A is chosen too small to 
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Fig. 6. Dynamics of Pi(t) for a constant interaction with Vij — 1, parameters: A'^ — 500, Se = 0.5, 
A = 0.00005. 



fulfill criteria ([5]). One again finds no complete decay into equilibrium here. However, one can 
see that Fermi's Golden Rule is fulfilled because of the good agreement at the very beginning. 
The deviation at some certain time nearly coincides with the Heisenberg time that is equivalent 
to the recurrence time of the autocorrelation function of V in our model with equidistant local 
energies {th ~ 6200). At this time the assumption of a completely decayed autocorrelation 
function necessarily becomes wrong and therefore the rate can no longer be considered as 
constant. Since the relaxation time of Pi{t) is increased by a smaller A, while the Heisenberg 
time is not, this finite size effect becomes important. 

Generally, the full evaluation of the second order Dyson term or even higher orders is nu- 
merically and analytically very extensive. Instead of that one can introduce a specific structural 
requirement for the interaction to estimate those contributions. In [5] Van Hove used a similar 
interaction structure to derive the possible occurrence of statistical relaxation. Basically, this 
would mean that V features this so called Van Hove-structure, if is dominated by its diagonal 
elements in some sense. If the interaction does not have Van Hove-structure, then statistical 
relaxation cannot emerge. If the interaction possesses Van Hove-structure, the dynamics are 
possibly but not inevitably statistical. So it can be regarded as a necessary criteria. One finds 
that the random interaction fulfills the Van Hove-structure, whereas the constant interaction 
does not. In this sense the resulting dynamics for both examples are explained correctly. In prin- 
ciple, the Van Hove-structure gives an estimation for the proper convergence of the respectively 
used perturbation expansion. It indicates that the higher orders should be small compared to 
the leading one, at least for times in the order of the relaxation time. 



3 Analytic derivation of the solution of the rate equation 

This Section is devoted to the derivation of the solution of the probability for the excitation 
to be on the left Pi{t) or right hand side P2{t) for the case where the interaction is modeled 
by complex Gaussian entries in V. We will not show the total derivation as some intermediary 
results are too long but rather try to explain the main steps. This is why some results are 
just stated. What we mean by solution is the following. If we keep the interaction V as a 
random matrix with its probability distribution instead of fixing it, then the time evolution 
operator, Q-^iHo+\v)t^ becomes a random evolution operator. Thus Pi{t) also contains this 
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random character. We can then calculate its average, E[Pi(i)] = / p(V)dVPi{t), where p{V) is 
the probability distribution over the random interaction. p{V) is the product of the probability 
distribution of the complex Gaussian elements which can be written in a more compact form 
as p{V) = ■^e^~'^''^^ \ where Z is the partition function. Wc will calculate this average in the 
limit where the number of energy levels per subunit (there are 2 subunits in the present case) N 
tends to oo and in the Van Hove limit. The Van Hove limit is given by lim^^o lim^^oo = T. 
We call T the macroscopic time. This limit represents some long time-weak coupling limit. We 
will actually not calculate E[Pi(t)] directly. We will calculate E[Pi{t) — P2{t)]- Since we know 
that Pi{t) + P2{t) = 1 always holds (may it be averaged or not), we can obtain E[Pi(f)] from 
it. 

The main idea of the proof and the below presented expansion is that, in the limits con- 
sidered, the interference effects vanish. For this graphs are introduced which represent pairs of 
histories of the state of the system. The vanishing of interference effects means that for two 
different pairs of histories the contribution is much smaller than for equal or similar pairs. In 
some sense the quantum feature of interference becomes negligible when the limits are taken, 
thus making it more classical. It is basically this fact that allows for an autonomous equation 
to exist. For a pair of histories, or a graph, it is the random interaction, after averaging, that 
will decide on the weight given to it. It will also decide which class of pairs of histories, or 
graphs, are preponderant, thus contributing the most. These will be called simple graphs. We 
start by writing down the observable we wish to consider and inserting the expansion of the 
time evolution operator. 

N 

P,{t) = (V^ole^^*^ |l,0(l,Z|e-^^*|Vo> (13) 

1=1 

Pi'(t) = (Vo|e'^*|l,0(l,«|e-^^*|Vo) (14) 

Every index I represents the dependency on Ei and \fi, I) is a basis ket for the state being in 
unit yLt on energy level Ei. 

Loosely speaking the idea is to expand the time evolution operator in powers of the interaction. 

We have then powers of the random matrix over which we can average. Some parts won't 
contribute in the limit N ^ oo and others won't contribute in the Van Hove limit. In the limit 
A'' — » oo we will keep our local spectrum bounded and so the energy level variables, Ei, will 
turn into continuous variables, E. Pi{t) will become then Pi{E,t). We expand the evolution 
operator as 



Pnit) 

Wc thus have 

pi{t) = Y.{t\r{ 

N 

n,m ;„,i^ = l p„ ,pj^ = l 

(18) 

We notice that if n + m is odd, (p^, l'„^\rj-,^\l, Z|_r„|p„, In) is proportional to an odd number 
of random variables with average zero. The average of an odd number of random variables 
centered around zero is zero and so all terms with odd n + m don't contribute to the average. 
For n + m to be even we need both to be even or both to be odd. The variables p'^ and p„ 



^(-zA)"r„(t) , 

n 

... dsQ... dsne-'"°''°V . . . Ve-'^^'-S *J j 



(15) 
(16) 



-ixr{UrL{m,mi\rn{tMo) , (17) 

2 

E ro{pLOMPnJn)i'^{-irX"+"'{pL,L\rl\l,l}{l,l\rn\Pn,ln) ■ 
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here stand for the unit 1 or 2. Since the interaction matrix makes an excitation hop from unit 
1 to unit 2 and vice versa, we notice that if n is even, then the product (1, l\rn\pn, In) is zero if 
Pn is equal to 2. Thus p„ has to be equal to 1. For n odd we will only have a contribution for 
Pn = 2. We make this explicit in formula (fT9| . 



N 

P^(t) =^A^("+™) E ro{lJ2rn)Mhl2n)^''''{-^r''{lJU^L\hl){lJ\^2n\l,l2n) 

N 

(2,Z2™+il^L+i|l,0(l,^l^2„+i|2,W> (19) 
-^iW + ^iW (20) 

Notice that the first contribution is related to the initial state of unit 1 and the second to the 
unit 2. Of course, we have a similar expression for ^2(^)1 

Pi{t)^Gl{t) + G[{t) . (21) 

When subtracting P2 from Pi we can group terms with the same initial data dependence that 
is FI with G[ , and with Gj. We then have 

E[Pi{t) - Piit)] = E[Fl{t) - G[{t)] + E[F^(t) - G^it)] . (22) 

We will thus compute one of these, the other being analogue. Inserting ([16]) in {l,l\rn\pn,ln) 
and identities after each interaction term we have 



{l,l\rn\Pn,ln) = j\ . . J^sq . . . ds„e-'^'^« e'^^'i . . . e-^^''.^"<5 I t - E j (23) 



X {~l)''{l,l\V\pi,h){pi,h\V\p2, h)--- {pn-lX-l\V\pn. In) 

For a shorter notation we define 



K-{t,{El^)^{-^) 



... dso...ds„e-'-^'"''e-*^'i''i ...e-^^'-^-c^ l^i-E^jj ,(24) 

i"(Oj},te}) = {1,1\V\P1M){PIM\V\P2M) ■ ■ ■ {Pn-uln-l\V\Pn,ln) ■ (25) 

All of the randomness is encoded in L"({Zj}, {pi}) and we know the {pi] variables are 
determined if n is odd or even. With this notation we have for the first term in Eq. ^9)1 

iV 

E V'0*(1, ^2™)V'0(1, Z2„)*'"(-*)'"(l, 4nl^2ljl, 0(1, Z|r2„|l, /2„) 

- E ^^5(1, Zi™)V'o(l, ;2„)^'"(i, {i?/.})^""(i, {i?/j})i'"({/,}, fe})Z2"({/;}, {P^}). (26) 

To average we must average over This term is a product of Gaussian complex variables. 

It only contributes if these variables correlate. This will introduce relations between the different 
indices {U, I'j, I}. According to Wigner's theorem we have the following formula to calculate the 
average over a product of Gaussian random variables. 

E[Xi...X2„] =E n (27) 
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Fig. 7. Example of graph representation for a pairing with n = 5 and m = 7. This graph represents 
the following contribution: 

m'7\VK){W,)] X E[{l'5\V\Q{l'2\V\l[)] X . . . E[ai|\/|0(Zi|\/l/2>] X E[(/2|\/jZ3>(/3iV|M]. 
All pj indices are omitted for lighter notation. 



Here Ctt is a permutation of the set of indices {1, . . . 2ti}. In the present case each Xj corresponds 
to matrix elements of the type {pi,h\V\p2, h)- The pairings induce a graph structure on L^L"^. 
Pictorially an example is given in Fig. 1. The average of one pairing is given in Eq. 



E[{l,h\V\2,h+,){2,l,\V\l,l,+,)] - ^Si^^i^^,Si^^,,i^ (28) 

Notice that the averaging enforces identities amongst {li,li+i,lj,lj+i} and hence amongst the 
energies {Ei., Ei^^-^, Ei. , i^ij+i} in Eq. The question is then how much graphs weight and 
what kinds of graphs are important. We focus on the weight they have. 

Over many pairings the delta functions are responsible for the graph structure. However, 
we notice the weight of every pair will always be the same. Thus we conclude the weight of a 
graph on L2«j^2m ^ gj^^j^ (w)"^"- E[i^"({?i}, {Pj})L^"'i{l'j}, {p'j})] is represented 
by a graph and enforces identities amongst the energies {Ei^, Ei' }. We can then rewrite p6p as 



:raphB C-n- (2n,2m) \ indcp.{/^ } 



/ (J \ n-f-m _ 

(^) roih l'2m)Mh l2n)Kl-{t, {Eu , E[;\)Kl^{t, {E^ , Ei.}) 



(29) 

Here K^'^it, {Ei^, E'iJ)Kl"'{t, {Ei^, Ey^}) stands for the function K'^'^K^^ from Eq.® but 
with the identities amongst the energy variables imposed. It can be shown that for C^{2n^ 2m) 
graphs the maximum of independent {li, 1'^} is n + m and that wc only have this for a certain 
class of graphs, which we will call simple graphs (S. graphs). For this class of graphs holds 
that the identity l2n = ^2m fulfilled. The fact that this identity is fulfilled implies that 

We would then have 



/ (J \ n+m / / 



n+m 



The Wj here represent the independent Ij and / dE is a short notation for the multiple integrals. 
The sums over independent variables turn into integrals over, what has become, continuous 
energy variables Ej . The limit as iV ^ cxo of the average of Eq. ([26|l then becomes 

nFi{Eo,t)] ^ lim E[F{(t)] (31) 
= $](A2)"+™ ( / dEP,{E2n. t = 0)if^(i, {E,])Kl"'it, {E,}) 

n,m s graphs C„(2n, 2m) 

Eq is the now continuous variable Ei. We now explain what type of graphs contribute in the 
way previously described. Among all graphs the ones of the type shown in Fig. [5] are simple 
graphs. Figure [5] represents graphs where n random variables from the left hand side are paired 
with one on the right hand side, and in between each pair of this left-right type of pairs we 
have a certain amount of pairings with their nearest neighboring random variables. That is. 
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(j + l)*'^ contraction 
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Fig. 8. Sketch of simple graphs 



in between the j*'' pairing from left to right and the {j + 1)*'' pairing from left to right we 
have kj pairings on the left and kj pairings on the right. It is very important to notice that 
the information of the graph is now coded in the set of number {n,kj,kj} which determine 

uniquely the graph. For such graphs the exponentials of the K"K"^ factors in Eq. ([26|l and Eq. 
PTjl in between the j*'' and (j + 1)*'* contraction have the form 

Cj depends on many variables ( Cj = Cj{Ej, Ej^q, E'j ^, , Sjq^ Tjq', sf, rf) ) which we didn't write 
explicitly. In order to calculate the contribution of a graph we need to integrate the product of 
all Cj functions over all of these variables respecting the delta functions of Eq. p4)) . It can be 
shown that in the Van Hove limit each integral over the functions depending on Ej^q, Sjq, E'j ^, 
and Tjqi can be replaced by a given constant. We have then 

CjiEj,slT^) = i-tp+^e-'''^ S'ir '10'^^ (i)'=Hie»^^ E'lV ^'©^^ (33) 

and 

Q{C^{n,{h,k'^}),Eo,t) = (^J dEP,{E2n,t = 0)Kl''{t,{E,})kl"'{t,{E,})^ (34) 



n / dE,P,{Er.,Q) n / dsj [] / dr/<5 \t -Y^sj ] 6 It -J^rn C,{E, 



j=i (=1 - " (=1 - " y J.J 

For the integration over the time variables we have the identity p5p . 

We have the same identity for the r variables. We see that with this identity we can make the 
integrations over the s\ variables and r/ variables easily because the exponential in Cj depends 

on Y^i^^ s\ and Y^iCx ■ 

'^J+l poo / \ '^Hl poo ( \ 

= (-i)'=.+i(i)*^-Hie-»^^-^i^0'=^e*^^"^-^^0'=^- (36) 

k'j\ 
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We then have by Eq. dM]), ([Ml) and jMl 

QiC^in,{k,,k'^}),Eo,t) 



(37) 



Fl{t) has become Fi{Eo, t) due to Ei — > i?o- Since our graphs are determined by the variables 
{n, kj, kj} and we have to sum over all possible graphs, this sum over all graphs becomes a sum 
over all possible values of {n, kj, k'^}. 

E[Fi(i?o,<)]= 5] (aA^f ^ (aA2)'='+^-^Q(a(n,{fc„fc;.}),£;o,t) (38) 

all ovon ft {ki,k'^} 

The fact that n is restricted to even numbers comes from Eq. (|20p. That means, Fi{Eo,t) is 
the contribution when we have an even number of interaction matrices on the left and on the 
right. This is guaranteed when n is even. Gi{Eo, t) turns out to be the sum over all odd n. By 
grouping together all terms which are to the power of kj and k'^ we can now sum over these. 

Y^{cjX^f^+'''.Q(C^(n,{h,k'j}),Eo,t) = (39) 

1=1-^ j=o ° 

In the first exponent — O is thus the imaginary part of 0. We have not derived this constant 
but it turns out to have the following imaginary part. 

0-0 = -i2n (40) 

Thus 

lim e-**"^'(^-^) = . (41) 

Van Hove 

We now turn to the first part of Eq. (P^ . By Eq. and Eq. ([M)) we have an expression 
for 'E[Fi{E,t)]. A similar one can be obtained for E[Gi{E,t)]. We group these together in the 
following equations. 

E[F,{E,t)-Gi{E,t)]= (^^')" E M2)'^-'+'=;Q(a(n,{fc„fc^.}),^,t) (42) 

all ovon ft {fci.fc^} 

- E (^^')'' E M2)'='+'=^Q(a(n,{fc„A:;}),£;,t) 

all odd n {ki,k'.} 

Up n 

^ ^^-.t.A^(e-e)-Q / dE,P,{Er.,0)l[ (43) 

all ft i=l •' j = 

i-aX'f ds.dr.S (t - ^ s,) ^ (t - E 

To calculate the integral part we introduce a change of variables, Eq. P^ . in the last line of 
Eq. O. 

«.=A^^,&.-^^ (44) 



3=0 



2 ' " 2 

+ o 
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In the Van Hove limit Eq. (j45|) turns into 

- rv — 1 
rpn ^ _^ 

(-'^)Vn2^'^(^^--^«) • (46) 

3=0 

Inserting Eq. in Eq. we obtain 

lim I dE(E[Fi{E,t)~Gi{E,t)]) 

Van Hove J 

= E ^"'"""^ n / 0)(-a)«— n 27r<5 - i?^) 

all ft j=r ■ j=0 

dEPi{E,0) . (47) 



all n 



In the same manner we can calculate / dE (E[^2(-E', t) — G2{E, t)]) in the Van Hove limit which 
will have the same form as Eq. ([47|l but will depend on P2{E, 0) instead of Pi{E, 0). Using these 
last results in Eq. (|22p we finally have 



lim lim (E[Pi(t)-P2(i)])=e-4^"^(Fi(0)-F2(0)) . (48) 

Van Hove Af^OO 

These give us the solutions to the rate equations with a rate of Ana. 



4 Conclusion 

The analytical derivation has shown that the dynamics of the occupation probabilities of single 
subunits are statistical for the average over the whole ensemble of possible interactions for 
N oo and in the Van Hove- limit. This means that if "most" members of the ensemble (a 
dense subset of the ensemble) give the same type of relaxation, then there is statistical relaxation 
which is then a typical feature for members of the ensemble. The numerical calculations for a 
finite size version of the system demonstrate that the occurrence of statistical relaxation actually 
depends on the structure of a concrete realization of the interaction. A random interaction 
reproduces the results for the average. Nevertheless, the example of the constant interaction 
indicates that there are exceptions as the dynamics are not statistical at all in this case. The 
fact that a fixed randomly chosen interaction reproduces the same type of relaxation as the 
average strongly supports that the statistical relaxation is a typical feature because it represents 
the majority of all possible interactions very well. 
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